Self— force gravitational waveforms for extreme and intermediate mass ratio inspirals. 
II: Importance of the second— order dissipative effect 
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We consider the importance of the second-order dissipative self force for gravitational wave 
dephasing for an extreme or intermediate mass ratio system moving along a quasi-circular 
Schwarzschild orbit. For the first-order self force we use the fully relativistic force in the Lorenz 
gauge. The second-order self force is modeled by its 3.5 post Newtonian counterpart. We evolve the 
system using the osculating orbits method, and obtain the gravitational waveforms, whose phase 
includes all the terms — within our approximation — that are independent of the system's mass 
ratio. The partial dephasing due to the second-order dissipative self force is substantially smaller 
than that of the first-order conservative self force, although they are both at the same order in the 
mass ratio. 

PACS numbers: 04.25.-g, 04.25. dg, 04.25. Nx, 04.70. Bw 

I. INTRODUCTION AND SUMMARY 

When a small compact object moves in the spacetime 
of a much more massive black hole it accelerates and 
emits gravitational waves. The cause of the acceleration 
of the smaller object fi may be viewed as the object's in- 
teraction with itself in the external field of the more mas- 
sive black hole M, namely the object's self force (hence- 
forth the SF). The SF causes a correction to geodesic 
motio in the mass ratio of the system e :— fi/M. The 
force on /i vanishes for geodesic motion. Its leading term 
in an expansion in e, that is the force at 0(e 2 ), has been 
studied in detail* 0, Q. For more detail, see the Intro- 
duction to Ref. [3| (hereafter Paper I) that will serve as 
an introduction to the present Paper too. 

It has been realized, however, that the dephasing 
of gravitational waveforms includes important contribu- 
tions also from the SF at 0(e 3 ) [4|]. We illustrate this 
point for quasi-circular equatorial orbits of (i around a 
Schwarzschild black hole M: The SF in this case can eas- 
ily be separated between its dissipative and conservative 
pieces; the dissipative piece is made of the temporal and 
azimuthal components of the self force, while the conser- 
vative piece is described by the radial component thereof. 
We schematically expand the phase of the gravitational 
waveforms in powers of e as follows: 



$ = e- 



$(0) +e$ (D +O(e 2 ) 



(1.1) 



For enumeration purposes it is useful to think in terms of the 
SF per unit mass. Then the leading order term is at O(e) and 
the next term at 0(e 2 ), which justifies the names first order and 
second order SF, respectively. As we refer here to the SF, not the 
SF per unit mass, the first order SF is at 0(e 2 ) and the second 
order SF is at 0(e 3 ). 



The dissipative SF at 0(e 2 ) causes the phase of the emit- 
ted gravitational wave to evolve at 0(e _1 ), that is it de- 
termines $(°). In essence, that is the same waveform 
as obtained when the orbital evolution is computed used 
balance arguments. The dephasing of the waveform at 
O(e ), given by $W, has two contributions: both the 
conservative piece of the SF at 0(e 2 ) and the dissipative 
piece of the SF at 0(e 3 ) contribute to the dephasing of 
the waveform at the same order, a dephasing which is 
independent of the system's mass ratio. 

It was shown in Paper I that the dephasing "jW at 
O(e ) of the gravitational waveform due to the first- 
order conservative piece of the SF at 0(e 2 ) is important 
for gravitational wave detection and astronomy. In this 
Paper we consider the other source of dephasing of the 
waveform at the same order, that is at 0(e ), specifically 
the second-order dissipative piece of the SF at 0(e 3 ). 

The SF we use in practice is the SF acting on a parti- 
cle moving for all eternity along a circular geodesic. As 
the particle does not move in actuality along an eternal 
circular orbit, but rather along a slowly decaying quasi- 
circular orbit, the actual SF differs from the SF obtained 
for circular geodesic by terms which are second order at 
0(e 3 ): the true SF is obtained from integration over the 
half infinite world line of the particle. Specifically, when 
one computes the leading order dissipative SF for the 
inspiralling world line, it would differ from its circular 
geodesic counterpart by terms at 0(e 3 ). We neglect such 
0(s 3 ) corrections to the SF: our computational method is 
the method of osculating orbits Q ; specifically, at each 
point along the orbit we assume that the orbit locally 
osculates to a geodesic, and it is the parameters of the 
geodesic that slowly evolve along the orbit. Therefore, 
the consistency of our approach allows us to consider the 
SF at each point along the orbit as that of an eternal cir- 
cular geodesic, so that in fact we include all O(e ) terms 



to the phase $ in our osculating orbits approach. 

Despite impressive progress in the understanding of the 
second-order SF (at 0(e 3 )) [f| 01: & i s as y et unavailable 
for us for practical use. We may however understand its 
significance qualitatively by mimicking it with its post 
Newtonian (PN) counterpart. (Our analysis can read- 
ily be adapted to the fully relativistic second order SF 
when the latter becomes available.) However, our con- 
clusions based on a hybrid model in which the first order 
SF is fully relativistic (specifically we use the Barack- 
Sago SF [8|]) and the second order is approximated by 
the post Newtonian expression will remain qualitatively 
unchanged when the latter is replaced by its fully rela- 
tivistic counterpart. We emphasize that our results are 
quantitatively accurate only in the far field, but close to 
the innermost stable circular orbit (ISCO) we lose accu- 
racy, because of our use of the post Newtonian expression 
for the second order SF. Specifically, the second order dis- 
sipative SF is known to 3.5 PN order, that is to 0(v 7 ) 
[9[. It is listed explicitly in Appendix [Al The first order 
dissipative self force at 3.5 PN overestimates the fully 
relativistic Barack-Sago self force by not more than 10% 
approaching the ISCO (and significantly less away from 
the ISCO). It is reasonable to expect that a similar ac- 
curacy is provided also by the second order dissipative 
force. Indeed, Isoyama et al [9] showed that the dephas- 
ing due to the 0(e 3 ) SF for quasi-circular extreme mass- 
ratio inspirals (EMRIs) may be well captured by the 3 
PN term. 

Unlike its scalar field counterpart, the gravitational SF 
is gauge dependent. Specifically, one can find a gauge in 
which the SF vanishes (a problem that would be equiv- 
alent to computing the SF driven orbital evolution), or 
gauges that satisfy certain algebraic or differential re- 
lations. The SF that we use in practice, namely the 
Barack-Sago SF [8J , is given in the Lorenz gauge. More- 
over, the orbit itself is also a gauge dependent quantity. 
Fixing the gauge to the Lorenz gauge, we focus attention 
on gauge invariant quantities ( "observables" ) , such as the 
waveforms and the plot of u , the temporal component of 
the four- velocity ("gravitational redshift," "helical Killing 
vector of the perturbed spacetime" ) as a function of the 
angular frequency fi. 

We first evolve the orbit and find the gravitational 
waveforms when the orbital evolution is driven by the 
first order SF (similarly to the analysis in Paper I), and 
then we add the contribution of the second order dissi- 
pative SF. We are mostly interested in the relative de- 
phasing of the two waveforms — that is, in the partial 
dephasing due to the second-order dissipative SF — and 
specifically in how large this dephasing is compared with 
the relative dephasing between the energy balance wave- 
form and the waveform obtained with the inclusion of 
the first order conservative SF (omitting the second or- 
der SF, as was done in Paper I) — the partial dephasing 
due to the first-order conservative SF. 

We find that in the parameter range we search, specif- 
ically the mass ratios 10 -5 -10 -2 the contribution to the 



dephasing because of the second order dissipative effect 
is much smaller than the contribution of the first order 
conservative effect, although they are both at the same 
order in that mass ratio (they are both independent of 
the mass ratio). For the orbits we describe in detail, 
specifically quasi-circular Schwarzschild orbits that start 
at tq = 8M and evolve down to the ISCO, the dephas- 
ing due to the second order dissipative effect amounts 
to 8% (with a negative sign) of the contribution of the 
first order conservative effect. In the gauge invariant plot 
showing u* as a function of the angular frequency $7, the 
effect of the second order dissipative self force is even 
smaller, and is manifested by a slower motion of the data 
point representing the system along the curve as the sys- 
tem is evolving. The smallness of the effect may suggest 
that it is unimportant for gravitational wave detection, 
although it may be importance for accurate gravitational 
wave astronomy, e.g., for precise parameter estimation. 
Of course the effect could be of the utmost importance 
for a chaotic system. Although chaos has been shown to 
exist for certain Kerr orbits with spinning particles, such 
orbits are not likely to be relevant for EMRI gravitational 
wave astronomy [10j. 

We comment that even though we include in our wave- 
forms the second-order dissipative SF in addition to the 
first-order conservative SF, so that all sources of de- 
phasing at order O(e ) are included, the resulting orbits 
and corresponding waveforms are not exactly the same 
as their self-consistent counterparts, except perhaps in 
the adiabatic limit. To obtain true self-consistent or- 
bits and corresponding waveforms one needs to simulta- 
neously obtain the instantaneous solution of the coupled 
SF integrated equations of motion and the perturbation 
equations. Self consistent evolutions have been obtained 
for scalar field SF 
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yet not for the gravitational SF. 
Although! self-consistent evolutions have many advan- 
tages, we use instead the fully relativistic SF that was 
obtained for exact circular orbits [8] within the oscu- 
lating orbits approach. In the adiabatic limit the error 
introduced by the latter method should become negligi- 
ble. Our waveforms may also serve to test self-consistent 
waveforms when the latter become available. 

The organization of this Paper is as follows: In Section 
HIl we describe the numerical method that allows us to 
evolve the orbit and to extract the corresponding gravita- 
tional waveforms with high accuracy and computational 
efficiency. In Section [TTT] we describe our simulations and 
find the waveforms, and finally in Section HVl we analyze 
the results and find the dephasing due to the effect of in- 
terest, and in Section [V] we describe the gauge invariant 
plot. 



II. THE NUMERICAL METHOD 

In our time-domain approach towards numerically 
modeling an EMRI system, there are two distinct compu- 
tations that need to be performed very accurately. First, 



is the computation of the orbital trajectory of the small 
object in the spacetime of the large black hole includ- 
ing the SF effects under consideration in this work. This 
has been performed using the osculating orbits method 
that we describe in Paper I. Since we showed in Paper I 
that this method is reliable, as the produced orbits are 
consistent with those obtained with the independent di- 
rect integration method, we use here the osculating orbits 
method exclusively. As shown in Paper I, our osculat- 
ing orbits code converges with 5 th order. Special care 
was taken to make certain of the accuracy of the trajec- 
tory throughout its rather long duration (on the order of 
10,000 cycles for e = 10 -5 ) in this work. 

Second, we use this trajectory generated in the previ- 
ous step as input to our time-domain Teukolsky EMRI 
code. This code, that numerically solves the inhomo- 
geneous Teukolsky equation with a particle-source, has 
been used for studying EMRI related problems for over 
a decade. Essentially, the Teukolsky EMRI code is lin- 
ear, hyperbolic (2+l)D PDE solver using a time-explicit, 
Lax-Wendroff, 2nd-order finite-difference numerical evo- 
lution scheme. The particle-source term on the right- 
hand-side of the equation requires some care for a nu- 
merical implementation. Details on the capabilities, ac- 
curacy, convergence and performance of the code appear 
at multiple places in the literature (see Ref. [12| for a 
recent review). 

In order to be able to efficiently perform the type of 
long evolutions (on the order of a billion time steps) be- 
ing presented in this work, we used an advanced high- 
performance version of the Teukolsky EMRI code, that 
is designed to execute on a large parallel supercomputer 
accelerated by general-purpose graphics processing units 
(GPUs). Details on the OpenCL-based parallel imple- 
mentation and careful measurements of gains in overall 
code performance can be found in Ref. [1 31 ) (and refer- 
ences therein) . This code also has an added compactihcd 
hyperboloidal layer to the outer portion of our computa- 
tional domain [12|]. This advancement allows us to map 
null infinity to the computational grid and thus, we are 
able to extract gravitational waveforms directly at null 
infinity. These technical advancements have helped to- 
wards improving the performance and accuracy of our 
time-domain Teukolsky EMRI code by several orders-of- 
magnitude over previous versions. This code currently 
performs at a level of accuracy on the scale of a hundreth 
of a percent while still maintaining a very high degree of 
computational efficiency [13| . 



III. THE PERTURBATIVE WAVEFORMS 
A. The perturbative expansion 

The perturbative approach that underlies Eq. (11.11) as- 
sumes that one may expand the force acting on the small 
object fi in powers of the mass ratio e, that the result- 
ing orbital evolution can be similarly expanded, and that 




r/M 



FIG. 1: The quantity e _1 /3 as a function of r/M. The long 
(red) curve is the curve based on the fully relativistic Barack- 
Sago results for //\>, and the shotter dotted (blue) curve rep- 
resents the first term in the post Newtonian expansion. 



these expansion properties carry over to the correspond- 
ing waveforms. For e — > this is surely the case, but 
any physical system has a finite value for e. In Paper 
I we considered 10~ 3 < e < 10 -2 for practical pur- 
poses: the number of orbits and the resulting evolution 
time and number of numerical steps needed grow with 
e . While certain astrophysical systems of interest in- 
deed have mass ratios in the range studied in Paper I, 
we would like to study also mass ratios at the order of 
e ~ 10~ 5 , and in addition also quantify the mass ratio 
range for which our perturbative expansion of the orbit 
and the corresponding waveform is appropriate. 

The first restriction, as discussed in Paper I, is that the 
adiabatic condition holds, or that 8 ^> e 1 ' 2 , where <5 := 
p — 6 — 2e, p is the orbit's semilatus rectum, and e is the 
orbital eccentricity. This condition does not guarantee, 
however, that the corrections to the orbital motion in 
powers of e drop sufficiently fast. 

The expansion of the SF in powers of e, /g F = 
fci) + £3 f(2) + ' ' ' further reveals the condition that 

•C 1. We show in Fig. [1] e~ x /i as a 
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function of r. Clearly this condition is met everywhere 
for our entire mass ratio range. At large distances the 
post Newtonian expression is 
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A perturbative expansion of SF driven quasi-circular 
Schwarzschild orbits was done in Ref. [4J. Specifically, 
the radial velocity was expanded in powers of e accord- 
ing to v = EVru + e 2 W(2) + •• •• Specifically, we define 
7 := £U(2)/v(i). The results of [J] imply that typical ex- 
pressions for the radial velocity terms are such that the 
condition that 7 <C 1 can be given by 
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FIG. 2: The quantity E~ lr y as a function of r/M. The thick 
(red) curve is the curve based on the fully relativistic Barack- 
Sago results, and the thin dotted (blue) curve represents the 
first three terms in the post Newtonian expansion. The insert 
shows their ratio. 



The scaling of fr ~ e 2 implies that 7 ~ e. We therefore 
plot e _1 7 as a function of r in Fig. [2j The far field 
behavior of 7 can be described by the post Newtonian 
expansion, specifically 
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The condition (|3.ip and Fig. [5] allow us, for each value of 
the mass ratio e, to determine the domain of validity of 
our perturbative expansion. 

We consider the perturbative expansion again below in 
Section UnCl 



B. The waveforms 

We produce three waveforms: first, the waveforms for 
the case that only first order dissipative terms are in- 
cluded, and all higher order terms are excluded ("the en- 
ergy balance waveform"; dephasing to 0(e _1 ) (hereafter 
WF-I); the waveforms for the case that in addition to 
the former the first order conservative term at O(e ) are 
included (WF-II); and the waveforms that in addition to 
the former also the second order dissipative term, also 
at 0(e°) (WF-III). The three different cases are summa- 
rized in Table HI 

Figures [3] and [4] show the h 22 mode of the waveform for 
particles of mass ratios 10 -3 and 10~ 5 , correspondingly, 
that start inspiraling from r$ — 8M down to the ISCO (in 
practice we let the particles get to 6.01M) as functions of 
the time. In all cases the three waveforms are in phase at 
t = 0, and the dephasing has not gone over 2ir radians. 

Figure \E\ shows the waveforms for the modes K[l x for 
the mass ratio e = 10~ 5 for the same orbit as shown in 
Fig. |U At t = all three waveforms are in phase for both 
polarization states. Other modes can be shown similarly. 



FIG. 3: The three waveforms for the mode h^ for the mass 
ratio £ = 10 -3 for a particle that starts at rrj = 8M and 
whose orbit decays down to the ISCO as function of et/M. 
The observer is positioned on the equatorial plane. The insert 
shows the last wavelength. We show WF-I with the dotted 
curve (black), WF-II with the red curve (dashed), and WF- 
III with the blue curve (solid). 




FIG. 4: Same as Fig. [3] for the mass ratio e = 10 
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FIG. 5: Same as Fig. [3] for the mass ratio e — 10 for the 
modes h\l x . 



TABLE I: Summary of the orbits and corresponding wave- 
forms. The columns describe the effect for the phase of the 
gravitation wave, and the rows are the three types of wave- 
forms that we consider. See the text for more detail. 





Dissipative 
at 0(e _1 ) 


Conservative 
at O(e ) 


Dissipative 
at O(e ) 


WF-I 
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- 


- 


WF-II 


+ 


+ 


- 


WF-III 


+ 


+ 


+ 



Figures EllH and [5] suggest that the dephasing between 
WF-III and WF-II is substantially smaller than the de- 
phasing between WF-II and WF-I. We study this con- 
clusion quantitatively below in Section [IVI 
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FIG. 6: The variable a as a function of distance. The solid 
(red) curve is the fully relativistic expression. The shorter 
dotted (blue) curve is the first post Newtonian term. 



C. The perturbative expansion revisited 

The smallness of the effect raises the question of 
whether the second order dissipative effect is negligible. 
In Ref. (j) this question was raised, and specifically a 
condition for neglecting the second order dissipative ef- 
fect was considered. This condition is given by 
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The scaling of the self force components with e implies 
that a is at O(e ), i.e., independent of the mass ratio. 
At great distances a can be given by its post Newtonian 
expression as 
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Figure [5] shows a as a function of r. Both at great 
distances and very close to the ISCO a is very small, 
specifically a — > 0~ for either r — > oo or r — > QM. The 
parameter \a\ is maximal around r ~ HAf , and is never 
greater than 0.043. As a is independent of the mass 
ratio, we conclude that justification for neglecting the 
second order dissipative effect can be found only for or- 
bits very close to the ISCO (or alternatively for orbits 
that are only very far from it). 



IV. THE DEPHASING 

We study next the dephasing between the three wave- 
forms that we presented above in Sec. IIIIBl We present 
in Fig. [7] the relative dephasings between the waveforms 
as functions of the time for the h 22 mode, and in Fig. [8] 
we plot the same for the h^ 4 mode. The dephasing of 
WF-I and WF-II was presented in Paper I. For ease of 



comparison, it is reproduced in Figs. 17181 We also show 
in Figs. EU the dephasing of WF-I and WF-III, and the 
dephasing of WF-I and the waveform by adding only 
the O(e ) contribution of the second-order dissipative 
effect to the first order dissipative case. We also show 
in Figs. I7|8I the ratio of the latter the the second case. 
Except for the expected numerical noise, this ratio is at 
a nearly constant level of 8% throughout the evolution. 
All the curves in Figs. 17181 are independent of e, and the 
curves for the various values of the mass ratio that we 
considered (in the range 10 -5 -10 -2 ) overlap each other 
in this figure. 
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FIG. 7: The dephasing at O(e ) between the different /ij 2 
waveforms as functions of et/M. When plotted against this 
variable, the figure is independent of the mass ratio e. The 
solid (blue) curve shows the dephasing between WF-I and 
WF-II, the dashed (red) curve shows the dephasing between 
WF-I and WF-II, and the dotted (black) curve shows the 
dephasing between WF-II and WF-III. The insert shows the 
ratio of the dotted to the solid curves. 

Notice that the contribution of the the second order 
dissipative effect is to reduce the overall dephasing of the 
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FIG. 8: Same as Fig. for the h\ A mode. Notice that the 
dephasing itself is greater than for the /ij 2 mode, as can be 
expected from the shorter wavelength of the hf 4 mode. 
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FIG. 9: Overlap integral for the full waveform WF-III with 
a sliding window taken from WF-II as a function of the win- 
dow's width, for the mass ratio e = 10~ 5 . The circles repre- 
sent the data points, and the solid curve is a quadratic fit. 



waveform at 0(e) with the waveform at 0(e 1 ). This 
result is indeed expected as the sign of fLs is opposite to 
the sign of ffa. 

The contribution of the dephasing at O(e ) is therefore 
dominated by the first order conservative effect. The con- 
tribution of the second order dissipative effect is to de- 
crease the total dephasing, by an amount of about 8% 
of the total O(e ) effect. We therefore suggest that for 
gravitational wave detection the second order dissipative 
effect is not very important, although it may be so for 
accurate gravitational wave astronomy, e.g., for accurate 
parameter estimation. Neglecting this effect would re- 
sult in an inconsistency in the determination of the mass 
ratio, specifically in increasing the error bars of the mea- 
surement. Figures [7] and [5] are our main result in the 
present Paper. 

The relatively small contribution to the total dephas- 
ing of the second order dissipative effect is further mani- 
fest in the overlap integral between the full waveform at 
O(e ), i.e., WF-III, and an integration window selected 
from the waveform that includes at O(e ) only the first 
order conservative effect (WF-II). We select the window 
from the end of the waveform, and find the overlap inte- 
gral as a function of the width of the integration window. 
In Fig. [UJwe plot the overlap integral as a function of the 
window width. Comparing Fig. [UJ with Fig. 12 of Paper 
I, we see that the contribution of the second order dissi- 
pative effect to the total overlap integral is at the order 
of 1(T 3 . 



V. THE GAUGE INVARIANT PLOT 

We showed in Paper I that when the component u* of 
the four-velocity is plotted against the angular frequency 
£1 for WF-I and WF-II, the two curves overlap. That is, 
when the orbit of the small mass fi is represented on this 



gauge invariant plot, the curve is unchanged (or at least, 
is changed only very little). However, we showed in Pa- 
per I that the speed by which the data point representing 
the system moves along the curve on the gauge invariant 
curve depends on whether one includes or excludes the 
first order conservative effect. Specifically, with the in- 
clusion of the first order conservative effect the system 
moves faster along the gauge invariant curve compared 
with the energy balance case. 

In Fig. [TUJ we present the gauge invariant curves for 
WF-II and for WF-III. On the scale of Fig. [TUJ the two 
curves completely overlap. However, there is a small dif- 
ference between them, at order 1CP 4 . However, a data 
point describing a system evolving according the WF-III 
moves slower along the curve than a data point describ- 
ing a system evolving according to WF-II. The effect of 
the second order dissipative term is to retard the orbital 
evolution. The effect, however, is so small, that it is hard 
to see how the gauge invariant plot could be used to dis- 
tinguish between templates that include or exclude the 
second order dissipative term. 
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FIG. 10: The gauge invariant plot for the mass ratio e — 
10 -5 , showing u* as a function of Q.M for two cases: for WF- 
II (solid blue curve, □) and for WF-III (dashed red curve, 
o). The two curves are indistinguishable on this figure, and 
are independent of e. We also show equally spaced (in time, 
with increments of 10, 000M) data points. The data points 
corresponding with WF-III advance slower along the curve 
than those corresponding with WF-II. The upper insert shows 
the different between the WF-II case (denoted here Mm) and 
the WF-III case (denoted here u%\) as a function of f2M, and 
the lower insert shows the same as a function of et/M (not 
gauge invariant!) for the same data points as in the main 
figure. 



Appendix A: The self force expressions 

We expand the self force as /g F = e 2 fK,+e 3 ffa +••-. 
The expressions we use for the self force are given by (see 
Paper I for more detail): 
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and the expansion parameters are given in Table [XXJ 
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